Fokker-Planck equation with memory: the cross over from balhstic to diffusive 
processes in many particle systems and incompressible media. 
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The unified description of diffusion processes that cross over from a balhstic behavior at short 
times to normal or anomalous diffusion (sub- or superdiffusion) at longer times is constructed on 
the basis of a non-Markovian generalization of the Fokker-Planck equation. The necessary non- 
Markovian kinetic coefficients are determined by the observable quantities (mean- and mean square 
displacements). Solutions of the non-Markovian equation describing diffusive processes in the phys- 
ical space are obtained. For long times these solutions agree with the predictions of continuous 
random walk theory; they are however much superior at shorter times when the effect of the ballis- 
tic behavior is crucial. 



I. INTRODUCTION. 

In spite of the fact that particle diffusion process are 
being studied for about two centuries three are still 
some subtle issues that need to be faced at the present 
time. In this paper we point out these subtle difficulties 
and offer ways to overcome them. 

The quantitative theory of diffusion processes begins 
in 1855 with the phenomenological solution of the diffu- 
sion problem by Pick. Pick employed an empirical defi- 
nition of the particle flux through a surface of a subvol- 
ume (Pick's first law) and the continuity equation which 
reffects the conservation of particles. This combination 
results in the diffusion equation (the Pick's second law) 
which defines the time evolution of the probability dis- 
tribution function (PDP) of particle concentration. The 
variance of this PDP grows in time according to 



Dt, 



(1) 



where D is the diffusion coefficient which in general de- 
pends on the particles and medium in which they diffuse. 
The equation solved by the PDP is the classical diffusion 
equation 



df{r,t) ^^d^f{r,t) 



dt dr^ 
The solution of this equation with the initial condition 



(2) 



f{r,t = 0) = 5{r), 
where S{x) is the Dirac delta- function, reads 
exp{-r^ / (ADt)) 



f{r,t) 



(3) 



(4) 



This function is normalized to unity J f{r,t)dr = 1 and 



its variance agrees with Eq. ([1]). 

The first subtle difficulty that one needs to pay atten- 
tion to is that this solution allows particles to arrive to 
arbitrary distances from the origin in finite, and even in- 
finitesimal time. Moreover, in many cases processes that 



exhibit the law ([1]) for times larger than some time tc 
possess a different behavior at short times, i.e. 



(5) 



This behavior is referred to as "ballistic" . Recently, due 
to the development of the measuring equipment, it has 
become possible to observe the ballistic motion and the 
transition to normal diffusion in accordance with Eq. ([T]) 
[9l-[T]|. Therefore, the correct theory for the PDP should 
satisfy the two asymptotic limits given by Eq. ([5]) and, 
in general, by Eq. ([ij simultaneously. 

The cure for both the infinite speed problem and for 
the ballistic regime in 1 dimension was proposed by Davy- 
dov [1^ . who introduced an explicit time interval tc of 
the mean free path. In his approach the diffusion prob- 
lem is reduced to a solution of the well known telegraph 
equation (see, e.g., [11] )■ In contrast to the parabolic na- 
ture of Pick's second law the telegraph equation is of a 
hyperbolic type, resulting in a propagation of the front 
of the PDP with finite velocity. We reiterate however 
that it was stressed in [l^] that the telegraph equation is 
correct only in the one-dimensional case. The equivalent 
treatment for 2 and 3 dimensional diffusion does not ex- 
ist in the literature. One of the main results of this paper 
(cf. Sect. |V| is a way to achieve the same in dimensions 
higher than 1. 

Another issue that needs to be discussed is that the 
classical diffusion law Eq. ([1]) is not universally obeyed 
in Nature. Since the classical work of Hurst on the 
stochastic discharge of reservoirs and rivers. Nature has 
offered us a large number of examples of diffusion pro- 
cesses which are 'anomalous' in the sense that an ob- 
servable X diffuses in time such that its variance time 
dependence is 



t > tc 



(6) 



where a ^ 1. We expect a ^ 1 when the diffusion 
steps are correlated, with persistence for a > 1 and anti- 
persistence for a < 1 [1]. This law is usually valid only 
at long times and for t <C ic we again may find ballis- 
tic behavior. All the issues discussed above for classical 
diffusion reappear in the context of anomalous diffusion. 
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both in 1 and in higher dimensions. The bulk of this 
paper will deal with establishing methods to achieve a 
consistent theory for the PDF that is valid for all times. 

As a preparation for more complex situations we re- 
view in Sect. |ll]the telegraph equation that regularizes 
all the issues raised above in 1 dimension. In Sect. IIIII 
we generalize the methodology of the telegraph equation 
to Fokker-Planck equations with arbitrary non-local ker- 
nels in time. In the same section we show how to relate 
this kernel to the observable mean square displacement 
of the diffusive process. In the next section IIVI we turn 
to the Langevin equation for a guidance how to compute 
the mean square displacement to determine uniquely the 
associated Fokker-Planck equation. In the next Section 
we generalize the methods used in this paper to the 3- 
dimensional case taking into account hydrodynamic in- 
teractions. The last section offers a summary and con- 
clusions. 



This solution looks as if the initial conditions propagate 
with the velocity c without changing the shape of the 
distribution. 

For an arbitrary value of tc the solution of the telegraph 
equation is given by three contributions (see, e.g., [l5| 
and the references cited therein) 

/(e,r) = e-Uir-O 

+ [//(e,r)+///(e,r)]e(T-0 ■ (11) 

Here the two dimensionless variables are the normalized 
time T = t/tc and the normalized length £, — r / ^Dtc- 
Q{x) is the Heaviside function. The variance of this dis- 
tribution reads 

(at = 2[r-(l-e-^)]. (12) 



II. REVIEW: THE TELEGRAPH EQUATION IN 
1 DIMENSION. 

In order to develop a diffusion model with possible fi- 
nite speed of propagation it is necessary to take into ac- 
count simultaneously the mean particle velocity c and 
the time length tc of the mean free path 12|. Then the 
diffusion process is defined by the following equation 



d'f{r,t) , dfir,t) d\fir,t) 



9t2 



dt 



D- 



dr 



(7) 



where the diffusion coefficient D = (?tc- 

In the limit ic — >■ the time evolution of the PDF is 
possible only if I? 7^ 0, e.g., the mean particle velocity 
c — >■ 00. In this case Eq. ([7]) is reduced to the classical 
diffusion equation ^ and the solution of this equation 
with the initial condition ([3]) is Eq. This solution is 
normalized to unity and its variance agrees with Eq. ([l} . 

The delta-function initial condition can be considered 
as a limit of a sequence of continues functions, for exam- 
ple 



5{x) 
5{x, a) 



lim 5{x^a) 



(8) 



Comparing with Eq. ([2]) we conclude that normal dif- 
fusion is an inverted process to the limit defined by Eq. 

, the initial delta- function spreads in space when time 
increases, keeping the center of mass at rest at the origin. 

In the limit tc ^ 00 Eq. ([7]) transforms to the wave 
equation 



9t2 



(9) 



Its solution with the initial condition defined by Eq. ^ 
reads 



f{r,t)=5{ct-r) 



(10) 



This form clearly crosses over from ballistic to normal 
diffusion as r increases. The interpretation of the three 
contributions is as follows: the first term corresponds 
to the damped ballistic propagation of the kind defined 
by Eq. (|10l) . The next two terms exist only inside the 
compact region that expands with the velocity c. The 
fimctions //(^,t) and fii{S,,T) are given by 



exp(-T/2) (I ^ 



(13) 



and 



exp{—T/2) 



h 



(14) 



where Iu{z) are the modified Bessel functions of the first 
kind. The asymptotic behavior of these functions is de- 
fined by I^{z)z^oo ~ exp{z) / ^/ 2tt z , therefore, in the long 
time limit r 3> 1 Eq. ([TT|) is reduced to the solution of 
the classical diffusion equation given by Eq. (g]). 

The various terms in the solution of the telegraph equa- 
tion, Eq. (ITlT) . are shown in Fig. [TJ Physically, the first 
term describes the part of the initial mass of diffusers at 
the origin that shoots out as a propagating solution which 
decreases exponentially in time, providing more and more 
mass to the diffusive part of the transport process. The 
function fi{£,, t)Q{t — ^) represents the spreading of dif- 
fusers whose mass peak at the origin, not participating in 
the propagation of the first function. Note that this solu- 
tion is zero at t = 0, but it is finite for any positive time. 
The third term ///(^, t)8(t— ^) is seen to first increase in 
weight and in amplitude and later to decrease. It stems 
from the existence of a finitely rapid front that does not 
allow the unphysical infinite speed of propagation that 
characterizes the solution Q . Of course all together the 
sum of the three contributions is normalized. 
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FIG. 1: The three contributions to the solution of the tele- 
graph equation PDF [Eq. (|lip ] for different time intervals r = 
0.5, 1, 2, 4, 8. propagating contribution is shown in the upper 
panel. The delta-function is graphically represented by nar- 
row Gaussians (see Eq. ([8])). The functions ,fi,ii{(,,T)Q{T — ^) 
(see Eq. (|13p and Eq. (|14|l are shown in the middle and lower 
panels. 



III. FROM THE TELEGRAPH EQUATION TO 
THE TIME-NONLOCAL FOKKER PLANCK 
EQUATION 

In order to motivate the use of the Fokker-Planck equa- 
tion with time non-local kernels we rederive now the tele- 
graph equation in the following way: let us begin with 
the continuity equation in the following form: 



dfir,t) dT{r,t) 



dt 



dr 



(15) 



Here T{r,t) is the par ticle flux. In its turn this flux sat- 
isfies the equation 



dT{r,t) r(r,t)-ro(r,t) 



dt 



to 



(16) 



where To(r,t) is the same as the first Pick's law 



dx 



The solution of Eq. (HH) is given by [13, [11 



tr 



exp 



t-t' 



Qj.2 



(17) 



(18) 



Substitution of this equation to Eq. ((T5|) yields the time 
nonlocal diffusion integro-differential equation 



df{r,t) D 



dt 



tc 



exp 



t-t' 



tc 



d'fir,t'] 

Qj.2 



At'. 



(19) 



After time differentiation this equation is reduced to the 
telegraph equation defined by Eq. ([T]). 

Observe now [iH] that Eq. (|19l) is a particular case of 
the time-nonlocal Fokker-Planck equation, 

where W{t) is the kernel responsible for the non-Pickian 
behavior of the diffusion process. The Laplace transform 
of the solution of this equation is given by [l^ 



f{r,s) 



W{s) 



P(r, 



W{s) 



(21) 



were the function P{r, s) is the Laplace transform of the 
solution of the auxiliary equation with the same initial 
condition 



(22) 



which is identical with Eq. The Laplace transform of 
the solution corresponding to the initial condition defined 
by Eq. (O follows from Eq. (g]) with = 1 



P{r, s) 



(23) 



where Kxjii^z) is the modified Bessel function of the third 
kind. Substitution of Eq. ^ to Eq. ^ yields the 
solution of Eq. (^0)) 



/(r-,s) 



2 r Ky2{r^s/W{s)) 



(24) 



The Laplace transform of an even moment of the dis- 
tribution is defined by 



j Rx, s) 



x^^da;, m=l,2. 



(25) 
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Substitution of Eq. ^ to Eq. ^ yield 

Wis)"" 



= 2mr{2m)- 



(26) 



For the moment of zero order it is necessary to take into 
account that \im„i-^Q iriT {2m) = 1/2, therefore, {r^)s — 
1/s, the normalization condition. 

For the variance (m = 1) Eq. ([26l) reads 



W{s) 



(27) 



The kernel of the time-nonlocal Fokker-Planck equation 
is defined by the mean square displacement and, therefore 
(see Eq. ([20])) the PDF is also defined by this quantity. 



IV. DETERMINING THE KERNEL FROM THE 
LANGEVIN EQUATION 

The conclusion of the last section is that in order to ob- 
tain the appropriate Fokker-Planck equation for a given 
process we need to determine the time dependence of the 
variance {r'^)t- One way to do so is to measure this mo- 
ment from experimental data. On the other hand, if this 
data is not available, or if one wants to derive this infor- 
mation from the physics of the problem, another starting 
point can be the generalized Langevin equation [l, [?[ . 

The standard Langenin equation is Newton's second 
law applied to a Brownian particle where the random 
force acting on a particle is taken into account. As was 
shown in [20|-|22| that the generalized Langevin equation 
is written in terms of a time non-local friction force: 



j{t-t')^dt' + R{t), (28) 



where R{t) is the random force component which is un- 
correlated with the velocity, having a zero mean {R{t)) = 
0. The autocorrelation of the random force is related to 
the kernel in Eq. (1^51) by the fluctuation-dissipation the- 
orem [H SI: 



(i?(W)) = ^7(^-0, 



(29) 



where m is the mass of a particle. For the kernel of a 
kind 7(t) = 7o^(i) Eq. (^5)) is transformed to the stan- 
dard Langevin equation. In the general case the Fourier 
transform of the random force correlation function is col- 
ored, for example cf. [23|V 

In dimensionless variables defined in the sectionllllwith 
the diffusion coefficient defined by 13 = ksTtc/m and 

oo 

1/tc = / 7(^)d< the Laplace transform of the mean square 



displacement follows from Eq. (see, e.g., [2J]) and 

reads 



s2(s -I- 7(s)) 



(30) 



where the dimensionless kernel satisfies in the time do- 

oo 

main J 7(T)dr = 1. One can see from this equation that 



if lim^^oo 7(s)/s the function ~ and the 
mean square displacement exhibit the ballistic behavior 

at short times ~ t^. The shape of the function 7(5) 
at small s is responsible for asymptotic behavior of the 
mean square displacement in the time domain. 

Substitution of Eq. ^ to Eq. (gT]) yields the kernel 



W{s) = 



1 



s + 7(5) 



(31) 



This equation settles the relation between the memory 
kernel of the time-nonlocal Fokker-Planck equation and 
the memory kernel of the Langevin equation. 

The solution of Eq. ([20]) given by Eq. with the 

memory kernel defined by Eq. (|31|) reads 



/(e,.) = #^[^+7(.)]^i4i^SM). 

' Vtt^^ (eVs(s + 7(s)))i/2 



Eq. (|32l) consists of two terms, 

/(^,s)=C[s + 7(s)]$(e's(s + 7(s))) 

where 

The Taylor expansion of Eq. is defined by 

(z - zoT. 



(32) 



(33) 



(34) 



n=0 



Z 



Z = Zo 



For the modified Bessel function of the third kind 

5" K,/2{y^) _ (-1)" i^i/2+„(yi) 



dz" 



-1/2 



2" 



-1/2+n 



(35) 



(36) 



Eq. ^ defines that z = ^^s{s+^{s)) ^ z^-f{^{s)/2)'^, 
where zq — ^^{s + (7(s)/2))^, therefore, Eq. ([551) reads 



2n 



5/n= 1 gi/2+n(e[g + 7(g)/2]) / g7(^) 

V^i^o"' K[s + 7(s)/2)]i/2+" \2V2^ 

(37) 

It is suitable to rewright Eq. ([55)1 in the following form 

= r[{s + j{s)/2)+j{s)/2]^{r's{s + ^{s))) 
= m,s) + fH{^,s). (38) 



Taking into account that 



TT e 
2 7i 



(39) 
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Eq. (I57|) can be written as 



e[s + 7(s)/2]) 



2n 



1 Xi/2+„(e[g + 7(^)/2]) / C7W \ ,4.^ 
Let linis_^oo 7(s) = 7o- Under this assumption the sum in 



Eq. (|40| asymptoticaUy tends go to zero and the function 
defined by Eq. (HIH) reads 



exp(-C[s + 7o/2]) 



^[5 + 70/2] 

and the function ///(f , s) from Eq. ([55)1 is given by 

ea;p(-^[s + 7o/2]) 



(41) 



7o- 



s + 70/2 

The inverse Laplace transform of Eq. (|42]) yields 



(42) 



(43) 



proceeding to limit lim^_^o t) — shows, that this 

contribution to the PDF has nothing with the initial con- 
dition and is responsible for diffusion of injecting particles 
during the transport process. 

The nonvanishing term in the function //(r, s) in the 
limit of large s is defined by 



exp(-e[s + 7o/2]). 



Its inverse Laplace transform yields 

//(e,0-e-^-<5(r-0. 



(45) 



Therefore, the part of the PDF under consideration 
which is responsible for the initial condition and the fur- 
ther impulse propagation is contained in the first sum- 
mand in Eq. (1551) . 

Now we can isolate from the PDF defined by Eq. ((32)) 
the part corresponding to the diffusion process of the 
particles which are lost by the propagating impulse 

hff{^,s) = /(e,s)-ea;p(-e[s + 7o/2]) 

= $(^,s)e"^^ (46) 

where the function $(r, s) is defined by Eq. (15^ and Eq. 
(l39l) and reads 



$(e,.s) = is + %s)) 



(47) 



The inverse Laplace transform of Eq. (H71) yields the 
diffusive part of the PDF in the time domain 



(48) 



The importance of this result is that the explicit Heavi- 
side function is taking upon itself the discontinuity in the 
function fdiffiS,, t). The function $(^, t — ^) in the time 
domain is a continuous function. If it does not have an 
analytical representation in closed form it can be evalu- 
ated numerically, for example using the direct integration 
method 25] . 

Summing together the results (|45p and (|48p in the time 
domain we get a general solution of the non-Markovian 
problem with a short-time behavior, in the form 



/(e,r) 



0- 



(49) 



One can see from this solution that the diffusion reparti- 
tion of the PDF occurs inside the spatial diffusion domain 
which increases in a deterministic way. The first term in 
Eq. corresponds to the propagating delta-function 

which is inherited from the initial conditions, and it lives 
decreasing in time at the edge of the ballistically expand- 
ing domain. 

In order to estimate the PDF in the long time limit it 
is necessary to suggest the form of the memory kernel of 
the Langevin equation at small s. The reasonable choice 
is given by 



7(s)s 



(50) 



One can see from Eq. ((30|) that under condition s"^ 



> 



(44) 

(i.e., a < 2) the inverse Laplace transform of Eq. (j30p 



yields the variance in the form of Eq. (|6]). Substitution 
of Eq. ([50]) to Eq. (|32p defines the Laplace transform of 
the PDF at small s 



-is" 



„l-ct/2 



(51) 



In the general case the inverse Laplace transform of Eq. 
(ICT) is given by the Fox functions (see, e.g. [26l|). For 
a = 1 in the time domain the PDF is defined by Eq. (jj]) . 
For the special case a = the inverse Laplace transform 
of Eq. ([5T|) is independent of time and reads 



(52) 



This result coincides with the PDF from [27(. 

Below we demonstrate the application of the time- 
nonlocal approach with explicit examples. 

A. Standard asymptotic diffusion 



For 7(5) = 1 in Eq. ([31]) the generalized Langevin 
equation is reduced to the standard one which predicts 
the mean square displacement in the form of Eq. (E 
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The memory kernel of the time-nonlocal Fokker-Planck 
equation reads 



1 



W{s) 

s -r -L 

In this case Eg. ipS)) reads 

/d.//(^,s) = $(C,u)e-«", 



(53) 



(54) 



where w = s + 1/2, therefore, the inverse Laplace trans- 
form is given by 



(55) 



The function $(<^,it) for the particular case under con- 
sideration follows from the definition given by Eq. (j47|) 



^{i,u) = u^{i,u)-i + -^{i,u), 



(56) 



where 



exp 



^v?- - 1/4 

From the initial value theorem follows that 
«'((^,t = 0) = /^m„^^^i^'(^,M) 

exv[- I i I [wVl- 1/(4^2) - u]^ 



(57) 



lim,, 



1, 



(58) 



Therefore, the inverse Laplace transform of Eq. ([56]) is 
given by 

$(^,t) = ^vI/(C,i) + i*(^,t) (59) 
The inverse Laplace transform of Eq. ([57]) is given by [l^ 



(60) 



where Io{z) is the modified Bessel function. 
Substitution of Eq. ^ to Eq. i^S^ yields 



(61) 



Substitution of this equation to Eq. adduces the 

solution to the result for the telegraph equation given by 
Eq. (HH). 

The diffusion part of this solution (the sum of graphics 
shown in the middle and the lower panels in Fig. [1} is 
displayed in Fig. [51 As was discuss above in the long time 




1 2 



FIG. 2: The diffusion part of the PDF defined by Eq. ([TTJ 
for time intervals r =.5, 1, 2, 4, 8. 



limit this part approaches the solution of the classical dif- 
fusion problem given by Eq. ^ In order to measure how 
fast is convergence it is convenient to estimate the time 
dependence of the kurtosis of the distribution defined by 



K = 



3. 



(62) 



For the Gauss distribution k = 0, for the pure ballistic 
propagation k = — 2, Moments are defined by Eq. ([26]) 
and for the kernel from Eq. (|53p the second moment 
is given by Eq. ([12]) Calculation of the fourth moment 
yields 



(e^), = 12[r2 + 6(1 - e-") - 2t(3 + e"")]. 



(63) 



FIG. 3: The time dependence of the kurtosis of the PDF 
corresponding to the memory kernel defined by Eq. (|53|l 

The kurtosis calculated with these moments is shown 
in Fig. [3] At short times the propagation is ballistic 
with following transition to the standard diffusion during 
a large time interval. 



B. Anomalous diffusion 

In order to join the short time ballistic and long time 
anomalous (see Eq. ([6]) and Eq. ([50])) the Laplace trans- 
form of the memory kernel of the Langevin equation can 
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FIG. 4: The diffusion part of the PDF defined by the mem- 
ory kernel of the Langevin equation (see Eq. (|64|l for different 
values of the parameter a. Superdiffusion (q = 3/2, upper 
panel), regular diffusion (a = 1, middle panel) and subdiffu- 
sion (a = 1/2, lower panel). Time intervals from the top to 
the bottom r =0.5, 1, 2, 4, 8. 



be defined by 



7(s) = (2-a)- 



(64) 



(1 + s)"-!' 
It follows from this equation that 

70 = 2 - a (65) 
and, therefore, the singular part of the PDF is given by 

fi{i,t) = lexp{~^-^^^T)5{\ i I ~r). (66) 

The diffusion part of the PDF was estimated numerically 
and is presented in Fig. 2] 

In a different way the kernel of the time-nonlocal 
Fokker-Planck equation can be defined by the mean 
square displacement given in the time domain. In [l^ 



■ c(=2 

■ c(=3/2 

■ ci=1 

■ ci=1/2 

■ n-O 



10 10 



t 



FIG. 5: The time dependence of the kurtosis of the PDF 
corresponding to the memory kernel defined by Eq. (|64p 



this quantity was proposed in a form interpolating the 
short time ballistic and long time anomalous behavior 



(67) 



where < a < 2 and Iq is the crossover characteristic 
time, a.t t <^to the law ([57)) describes the ballistic regime 
and at t ^ to the fractional diffusion. 

Introduce now dimensionless variables {£,'^)t — 
{Ar'^)t/{2Dat^) and r = t/to. With these variables the 
last equation reads 



(1 



(68) 



The Laplace transform of Eq. (|68l) is given by 



ie), = {--l)-+{{a-l){--2) + s) ^r{a - 1, .), 

(69) 

where r(a, s) is the incomplete gamma function. At large 
s Eq. can be written as 



270 

^2 



(70) 



Substitution of Eq. (|69p to Eq. (|70p yields the estimation 
of the asymptotic value of the memory kernel 70 = 3(2 — 
a) which defines the time evolution of the singular part 
of the PDF. The results of numerical calculations of the 
diffusive part are shown in Fig. |6] 

The time evolution of the kurtosises of these distribu- 
tions is displayed in Fig. d?]). 

The reader should appreciate the tremendous role of 
memory. The ballistic part which is represented by the 
advancing and reducing delta-function sends backwards 
the probability that it loses due to the exponential decay 
seen in Fig. [U upper panel. Thus the diffusive part of 
the PDF is replenished with time and asymptotically ap- 
proaches the Markovian distribution. The speed of con- 
vergence and the evolution of the PDF shape depends 
on the interim behavior of the mean square displace- 
ment (or, that is the same, of the memory kernel of the 
Langevin equation). 
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V. HYDRODYNAMIC THEORY OF THE 
TIME-NONLOCAL DIFFUSION 

A. Memory kernel 

The friction force acting upon a spherical particle has 
been derived in [2^ [s^l (the modern discussion one can 
be found in [sT'l ) and is given by 



Ffrit) 




(71) 



where Ffr{t) is the friction force, R is the particle radius, 
and p and 77 are the density and viscosity of the solvent. 
Substitution of this force in the Langevin equation yields 
the fractional equation and the Laplace transform of the 
memory kernel can be found in a straightforward way 
[33 |. For simplicity we will use instead of this approach 
results obtained [16] for the mean square displacement. 

A system with hydrodynamic memory has two charac- 
teristic times 



9 Ms 



(72) 



where Mg is the mass of a solvent particle, S = GirriR is 
the Stokes friction coefficient and 



M + Ms/2 



FIG. 6: The diffusion part of the PDF corresponding to the 
mean square displacement given by Eq. (|67|) for different 
values of the parameter a. Superdiffusion (q = 3/2, upper 
panel), regular diffusion {a — 1, middle panel) and subdiffu- 
sion {a = 1/2, lower panel). Time intervals from the top to 
the bottom t =0.5, 1, 2, 4, 8. 



1 M \ 
- 1 + 2 \ tr 

9 Ms 



(73) 



One can see from Eq. (f73|) that tf > 0, therefore, it is 
reasonable to introduce dimensionless variables r — t/rp 
and P = Tfijrp. Then the projection of the mean square 
displacement on one of the axis reads 




FIG. 7: The time dependence of the kurtosis of the PDF 
corresponding to the mean square displacement defined by 
Eq. dSZl 



(Ax^(r)) =2DTF<T-2W^r + (/?-l) 



e"i'^er/c(aiy'r) 



ai — a-i 



(74) 



where D is the diffusion coefficient, erfc{z) is the com- 
plimentary error function and coefficients ai and 0:2 are 
given by 



ai,2 



(75) 



The dependence of the coefficients ai and 02 on the pa- 
rameter [3 is shown in Fig. [5] 
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FIG. 8: Dependence of coefficients cn,2 on the parameter /3. 



At short times r ^ 1 the mean square displacement is 
given by the Taylor expansion of Eq. (|74)) 



105 V TT 



16 I3t 



1 

24 



/?(/?- 3) + 1 



(76) 



One can see that the mean square displacement defined 
by Eq. (l76t corresponds to the ballistic motion. 

In the long time limit the mean square displacement is 
given by 



{Ax^{t)) = 2Dtft{ 1 ~ 2 



13-2 (/3-l)(/3-3) 




2V^T 



(77) 



Eq. ([77)1 agrees with standard diffusive. 

On can see in Fig. ([5]) that there are a few particular 
points. 

1. /3 — 0, {Ms — 0). Substitution of this value of the 
parameter j3 to Eq. ((75|) yields ai.2 — ±i and Eq. ([74)) 
is reduced to Eq. ([T2l) . 

2. /3 = 3 (M = Ms). In this case ai,2 = (^(3) T i)/2 
and the mean square displacement is defined by the error 
function of complex argument. 

3. /3 = 4 {M/Ms =5/8). Substitution of /? = 4 into 
Eq. dZll) yields ai(2) = 1. It follows from Eq. that 



lim {x (r)) 

Qi— >a2— 1 




(78) 



4. /3 = 9 (Af = 0). In this case the real parameters 
ai,2 = (3t\/5)/2. 



The mean square displacements for different values of 
the parameter (3 listed above are shown in Figl2] The hy- 
drodynamic interaction results in time delay of the mean 
square displacement, the larger value of the parameter /3 
the slow convergence to the asymptotic limit. 

It is necessary to note that the predictions of Eq^ ([74)) 
are in a good agreement with experimental data [9| 

The Laplace transform of the mean square displace- 
ment given by Eq. ()74p reads 



(Aeis)) 



s2(s + i + v;8i)' 



(79) 



therefore, the memory kernel of the Langevin equation is 
given by 



7(s) = 1 + ^/Ps. 



(80) 

Substitution of Eq. ([80]) to Eq. ^ yields the kernel of 
the time-nonlocal Fokker-Planck equation with due re- 
gard for hydrodynamics effects 

W{s) = , ^ ^ . (81) 



s + 1 




FIG. 9: Time dependence of the mean square displace- 
ment for different values of the parameter fi {{^{tY') = 
{x\r))/{2Drp)). 



B. Probability density function 

Solutions of the telegraph equation in higher than one 
dimension were considered in [33i] and it was shown that 
the result is unphysical, in some regions the PDF be- 
comes negative. Nevertheless, a reasonable PDF in ex- 
plicit form was obtained for two-dimensional systems in 
the frame of a persistent random walk model [3J| . There- 
fore, the possibility to describe high dimension diffu- 
sion process with finite speed of propagation with time- 
nonlocal approach is an open question. 

In the general case Eq. ([20)) reads in dimensionless 
units 



dr 



W^(T-T')Arf/(^,T')dr' 



(82) 
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where d is the space dimension and the Laplace operator 
is given by 



1 


dC 

d- 


-1^ 

1 d 


de 







(83) 



For pure baUistic propagation the kernel in Eq. ((82|) is 
defined by W{t) = 1. The time derivative of Eq. 
with this kernel yields d-dimensional wave equation 



Ad/(C,T) =0. 



(84) 



It is clear that the propagation of the ballistic impulse is 
specified by the function 



6{t-0 



d-l 



(85) 



where the normalization condition is defined by 

00 

/ /(Ci 1- Nevertheless, this function is the 



solution of Eq. (|M| the one-dimensional case only; for 
higher dimensions the solution differs from the function 
defined by Eq. (j85p . turning negative in some regions. 
This is the source of the failure of the telegraph equation 
for d > 1. Substitution of Eq. to Eq. (|M1) yields 



V5t2 



Ad 



Sir - 



1 



d_ 



d-2\S{T^^) 



(86) 

Therefore, we can suggest that the auxiliary equation 
which corresponds to Eq. ()86|) is given by 



(87) 



5/(e,r) _ d^f{^,T) , d-ldf{^,T) 



dr de 

{d~l){d-2) 



e 



fier). 



The Laplace transform of Eq. (I87t reads 

idfies) , f{d-i){d-2) 



oe 



e 



-s]fiC,s). 



This equation can be reduced to the Bessel equation [35] 
and its normalized solution is given by 



2-rf -^l/2(C7s) 



(89) 



The inverse Laplace transform of Eq. 
in the time domain 



yields the PDF 



/(e,r) 



-«V(4r) 



(90) 



The PDF given by this equation differs from that of stan- 
dard diffusion. Nevertheless, the second moment of this 





FIG. 10: The PDF defined by Eq. 1^ ior t ^ 0.5,1,2,4. 
From the top to the bottom P — 3,4, 9. 



functions is in agreement with Eq. ([T]) and the kurtosis 
of this distribution is given by k = 3. 

We suggest that the time-nonlocal Fokker-Planck 
equation is given by 



dr 



WiT-T')Vfi^,T')dT' 



(91) 



where the operator V is defined by 

■n - ^ ^ 0^-1 9 (rf-l)(d-2) 

ae ^ e 

The formal solution of Eq. ([OT]) follows from Eq. 



(92) 
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and Eq. (|5T|) 



the time dependence of the mean square displacement of 



-7(s) 



eVKs + 7(s)) 



1/2 



(93) 

This equation is similar to Eq. (|32| to within the coef- 
ficient which depends on ^. Therefore the nonvanishing 
part of the PDF for d = 3 and the memory kernel defined 
by Eq. ^ is given (see Eq. (gO]) and Eq. g!])) by 



1 



//(^,s) = ^exp 



r{s + 



1 



(94) 



For /? = the inverse Laplace transform of this equation 
yields the damped impulse propagation corresponding to 
Eq. ([55]) . In a general case the inverse Laplace transform 
of Eq. ([94]) reads 



4^ 



/3eV(16[r-e]) 



(95) 



In the limit r — >■ the function defined by Eq. ((95|) 
tends to three-dimensional delta-function, the initial con- 
dition of the problem under consideration. In infinitesi- 
mal time the delta-function due to hydrodynamic inter- 
actions becomes smooth and the PDF is free from the 
singular contribution. Therefore it is not necessary to 
pull out the ballistic contribution from the full PDF and 
the solution of the problem is given by 



/(C,t) = $(e,r-oe(T-o, 



(96) 



where the inverse Laplace transform of the continuous 
function 




X ea:p(- I i I •s[W 1 



(97) 



This can be estimated numerically. Results of calcula- 
tions for different values of the parameter /3 are shown in 
Fig. [TUl Note that a similar ballistic peak of the PDF 
was observed in simulations of test particle transported 
in pseudoturbulent fields js^. The kurtosis of the PDF 
for different values of the parameter /? are compared in 

Fig. [n 



VI. CONCLUSION 

The strategy followed in this paper is to construct a 
time-nonlocal Fokker-Planck equation which reproduces 




FIG. 11: The time dependence of the kurtosis of the PDF 
corresponding to the memory kernel defined by Eq. (|81|l 



an underlying process throughout the time domain, it 
should be stressed that the same mean square displace- 
ment may correspond, in general, to different models for 
the time-dependence of the PDF. Thus the predictions 
of the time non-local Fokker-Planck approach must be 
compared to experimental data for the PDF to assess its 
scope of validity and the quality of the approximation. 
The advantage of the proposed model is that it can deal 
with diffusion processes that cross-over from a ballistic 
to a fractional behavior when time increases from short 
to long values respectively. 

The general one-dimensional solution (1491) demon- 
strates the effect of the temporal memory in the form 
of a partition of the probability distribution function in- 
side a growing spatial domain which increases in a de- 
terministic way. The approach provides a solution that 
is exists at all times, and in particular is free from the 
instantaneous action puzzle. An extension of the em- 
ployed approach for higher spatial dimensions is used to 
study the implications of hydrodynamic interactions on 
the shape of the PDF. It is shown that singular ballis- 
tic contribution to the PDF is smoothed out during the 
propagation. The expansion given by Eq. ([57]) splits the 
Laplace transform of the PDF into two parts. One is the 
ballistic part which is the solution of the inhomogeneous 
differential equation in spatial variable with the initial 
condition of the problem. The second is the diffusion 
part, which solved the homogeneous equation and is zero 
at the initial time. In its turn the diffusion part is divided 
again into a few terms which correspond, in general, to 
different auxiliary equations. Therefore, each diffusion 
term in s— space can be multiplied by a coefficient C(s) 
so that the normalization condition and the mean square 
displacement are preserved. Also here there is a freedom 
of modeling, and this freedom was not exhausted in this 
paper. Further analysis and comparison to experimental 
systems are necessary to nail down the various options 
of modeling different physical realizations. This research 
will the subject of future publications. 
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